/**
 * This file is part of Graph-Tool.
 * 
 * Graph-Tool is free software: you can redistribute it and/or modify it under
 * the terms of the GNU General Public License as published by the Free Software
 * Foundation, either version 3 of the License, or (at your option) any later
 * version.
 * 
 * Graph-Tool is distributed in the hope that it will be useful, but WITHOUT ANY
 * WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
 * A PARTICULAR PURPOSE. See the GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License along with
 * Graph-Tool. If not, see <http://www.gnu.org/licenses/>.
 */
package graphtool.math.matrices;

/**
 * @author a.lunkeit - alunkeit(at)gmail(dot)com
 * 
 * This class represents a Matrix of numeric values. It provides basic functionality
 * like transposition, determinants etc. In this project, NMatrix is used to set up
 * the required functionality for working with graphs.
 */
public class NMatrix
{
  /**
   * The values building the matrix
   */
  double[][] _values;

  /**
   * The dimensions of the matrix
   */
  Dimension _dimension = new Dimension();

  /**
   * Creates a new matrix object from a given 2-dimensional double array
   * 
   * @param values
   */
  public NMatrix( double[][] values ) throws MatrixPropertyException
  {
    if( null == values )
      throw new MatrixPropertyException( "null array is not allowed for this object type" );
    
    this._values = values.clone();

    this._dimension.rows = values.length;

    if( this._dimension.rows > 0 )
      this._dimension.columns = values[0].length;
  }

  /**
   * Creates new matrix with all values initialized by 0.0
   * 
   * @param rows
   * @param cols
   */
  public NMatrix( int rows, int cols )
  {
    this._values = new double[rows][cols];
    this._dimension.rows = rows;
    this._dimension.columns = cols;
  }
  
  /**
   * Returns the count of rows in the matrix
   * 
   * @return
   */
  public int rows()
  {
    return _dimension.rows;
  }
  
  /**
   * Returns the count of columns in the matrix
   */
  public int columns()
  {
    return _dimension.columns;
  }

  /**
   * Set a value in the matrix
   * 
   * @param column
   * @param row
   * @param value
   * @throws IndexException
   */
  public void setValueAt( int row, int column, double value ) throws IndexException
  {
    try
    {
      if( column <= _dimension.columns && row <= _dimension.rows )
        this._values[row][column] = value;
      else
        throw new IndexException( _dimension, column, row );
    }
    catch( ArrayIndexOutOfBoundsException e )
    {
      throw new IndexException( _dimension, column, row );
    }
  }

  /**
   * Returns a value in the matrix
   * 
   * @param column
   * @param row
   * @return
   * @throws IndexException
   */
  public double getValueAt( int row, int column ) throws IndexException
  {
    if( column <= _dimension.columns && row <= _dimension.rows )
      return this._values[row][column];
    else
      throw new IndexException( _dimension, column, row );
  }

  /**
   * Returns the dimension of the matrix
   * 
   * @return
   */
  public Dimension getDimension()
  {
    return this._dimension;
  }

  /**
   * Returns a string representation of the matrix
   */
  public String toString()
  {
    StringBuffer buffer = new StringBuffer();

    for( int r = 0; r < _dimension.rows; r++ )
    {
      for( int c = 0; c < _dimension.columns; c++ )
      {
        buffer.append( String.format( "%+.2f ", _values[r][c] ));
      }
      buffer.append( "\r\n" );
    }

    return buffer.toString();
  }

  /**
   * Returns the transposed matrix
   * 
   * @return
   */
  public NMatrix transposed() throws MatrixPropertyException
  {
    double transposed[][] = new double[this._dimension.columns][this._dimension.rows];

    for( int cC = 0; cC < _dimension.columns; cC++ )
    {
      for( int rC = 0; rC < _dimension.rows; rC++ )
      {
        transposed[ cC ][ rC ] = _values[ rC ][ cC ];
      }
    }

    return new NMatrix( transposed );
  }

  /**
   * Returns true if the matrix is squarish.
   * 
   * @return
   */
  public boolean isSquarish()
  {
    return _dimension.columns == _dimension.rows;
  }

  /**
   * Returns true if the matrix is symmetric. In case that a matrix is
   * symmetric, the original matrix and the transposed matrix are squarish and
   * contain the same values. So this function does the transposing and the
   * comparison of both objects.
   * 
   * @return true, in case that both are the same
   * @throws MatrixPropertyException 
   */
  public boolean isSymmetric() throws MatrixPropertyException
  {
    if( !isSquarish() )
      return false;

    NMatrix transposed = this.transposed();

    return this.equals( transposed );
  }

  /**
   * Checks if both matrices are the same (this-object and m)
   * 
   * @param m
   *          - The matrice to be compared with
   * @return - true in case of successful check
   */
  public boolean equals( NMatrix m )
  {
    if( !this._dimension.equals( m._dimension ) )
      return false;

    // check element wise
    for( int columns = 0; columns < _dimension.columns; columns++ )
    {
      for( int rows = 0; rows < _dimension.rows; rows++ )
      {
        if( _values[rows][columns] != m._values[rows][columns] )
          return false;
      }
    }

    return true;
  }

  /**
   * Returns a clone of the matrix object
   */
  public NMatrix clone()
  {
    try
    {
      return new NMatrix( this._values );
    }
    catch( MatrixPropertyException e )
    {
      return null;
    }
  }

  /**
   * Adds two matrices. Returns a new object with the result
   */
  public NMatrix add( NMatrix matrix2 ) throws MatrixPropertyException,
      IndexException
  {
    
    if( !this._dimension.equals( matrix2.getDimension() ) )
      throw new MatrixPropertyException( "Two matrices should be the same dimension." );
    
    NMatrix sumMatrix = new NMatrix( _dimension.rows, _dimension.columns );

    for( int i = 0; i < _dimension.rows; i++ )
    {
      for( int j = 0; j < _dimension.columns; j++ )
      {
        sumMatrix.setValueAt( i, j, _values[i][j] + matrix2.getValueAt( i, j ) );
      }
    }
    return sumMatrix;
  }

  /**
   * Matrix multiplication using Falk Scheme: this * M
   * 
   * @param matrix2
   *          - The Matrix to be multiplied with this.
   * @return The result of the multiplication.
   * 
   * @throws IndexException
   * @throws MatrixPropertyException 
   */
  public NMatrix multiply( NMatrix matrix2 ) throws IndexException, MatrixPropertyException
  {
    double io[][] = new double[_dimension.rows][matrix2.getDimension().columns];
    NMatrix mr = new NMatrix( io );

    // Counter over the column count of Matrix m2
    for( int columnIndexM2 = 0; columnIndexM2 < matrix2.getDimension().columns; columnIndexM2++ )
    {
      // iterate through the first row and multiply the values with
      // the values of the first column of m2
      for( int rowIndexM1 = 0; rowIndexM1 < _dimension.rows; rowIndexM1++ )
      {
        double sum = new Double( 0 );

        for( int columnIndexM1 = 0; columnIndexM1 < _dimension.columns; columnIndexM1++ )
        {
          double m1Value = getValueAt( rowIndexM1, columnIndexM1 );
          double m2Value = matrix2.getValueAt( columnIndexM1, columnIndexM2 );

          sum += m1Value * m2Value;

          mr.setValueAt( rowIndexM1, columnIndexM2, sum );
        }
      }
    }

    return mr;
  }

  /**
   * Multiply all values of the matrix with a scalar value
   * 
   * @param scalar
   * @return A new matrix containing the result of the operation
   * @throws MatrixPropertyException 
   */
  public NMatrix multiplyScalar( double scalar ) throws MatrixPropertyException
  {
    double values[][] = new double[_dimension.rows ][_dimension.columns ];

    for( int i = 0; i < _dimension.rows; i++ )
    {
      for( int j = 0; j < _dimension.columns; j++ )
      {
        values[i][j] = _values[i][j] * scalar;
      }
    }

    return new NMatrix( values );
  }
  
  /**
   * Subtract a matrix, return the resulting matrix
   * 
   * @param matrix2
   * @return
   * @throws MatrixPropertyException
   * @throws IndexException
   */
  public NMatrix subtract( NMatrix matrix2) throws  MatrixPropertyException, IndexException 
  {
    return add(matrix2.multiplyScalar(-1));
  }
  
  /**
   * Expand the matrix by the given number of rows and columns. Returns a new object of type
   * NMatrix.
   * 
   * @param cols
   * @throws MatrixPropertyException 
   */
  public NMatrix symXpand( int cols ) throws MatrixPropertyException
  {
    if( cols <= 0 )
      return this.clone();
    
    double values[][] = new double[ _dimension.rows + cols ][ _dimension.columns + cols ];
    
    for (int i = 0; i < _values.length; i++) 
    {
      System.arraycopy( _values[i], 0, values[i], 0, _values[i].length );
    }
    
    return new NMatrix( values );
  }
  
  /**
   * Expand the matrix by the given number of rows. New elements are initialized by 0
   * 
   * @param rows Count of rows
   * @return The new extended Matrix object
   * @throws MatrixPropertyException 
   */
  public NMatrix xpandRows( int rows ) throws MatrixPropertyException
  {
    if( rows <= 0 )
      return this.clone();
    
    double values[][] = new double[ _dimension.rows + rows ][ _dimension.columns ];
    
    for (int i = 0; i < _values.length; i++) 
    {
      System.arraycopy( _values[i], 0, values[i], 0, _values[i].length );
    }
    
    return new NMatrix( values );
  }
  
  /**
   * Expand the matrix with the given number of columns. New elements are initialized by 0
   * 
   * @param cols The number of columns used for expansion
   * @return The new extended Matrix object
   * @throws MatrixPropertyException 
   */
  public NMatrix xpandCols( int cols ) throws MatrixPropertyException
  {
    if( cols <= 0 )
      return this.clone();
    
    double values[][] = new double[ _dimension.rows ][ _dimension.columns + cols ];
    
    for (int i = 0; i < _values.length; i++) 
    {
      System.arraycopy( _values[i], 0, values[i], 0, _values[i].length );
    }
    
    return new NMatrix( values );
  }
  
  /**
   * Reduces the matrix by the given count of columns and rows
   * 
   * @param cnt The count of rows and values to be used for symmetric reduction. If the
   * given value is negative, the matrix will be expanded
   * 
   * @return A new Matrix object. The call has no effect on <code>this</code>
   * @throws MatrixPropertyException 
   */
  public NMatrix symReduce( int cnt ) throws MatrixPropertyException
  {
    if( cnt == 0 )
      return this.clone();
    
    if( cnt < 0 )
      return symXpand( Math.abs( cnt ));
    
    double values[][] = new double[ _dimension.rows - cnt ][ _dimension.columns - cnt ];
    
    for( int i = cnt; i < _values.length; i++ )
    {
      System.arraycopy( _values[ i ], cnt, values[ i - cnt ], 0, _values[i].length - cnt);
    }
    
    return new NMatrix( values );
  }
  
  /**
   * 
   * @param iRow
   * @param iCol
   * @return
   * @throws MatrixPropertyException 
   */
  public NMatrix remove( int iRow, int iCol ) throws MatrixPropertyException
  {
    double d[][] = new double[ _dimension.rows -1 ][ _dimension.columns -1 ];
    
    int rI = 0; // counting index of reduced rows
    int rC = 0; // counting index of reduced columns

    for( int i = 0; i < _dimension.rows; i++ )
    {
      if( i == iRow )
        continue;
      
      for( int j = 0; j < _dimension.columns; j++ )
      {
        if( j == iCol )
          continue;
        
        d[ rI ][ rC ] = _values[ i ][ j ];
        rC++;
      }
      
      rI++;
      rC = 0;
    }
    
    return new NMatrix( d );
  }
  
  /**
   * Computes the determinant of a given Matrix.
   * 
   * @return
   * @throws NoSquareException
   * @throws IndexException
   * @throws MatrixPropertyException 
   */
  public double determinant() throws 
      IndexException, MatrixPropertyException
  {

    if( !isSquarish())
      throw new MatrixPropertyException( "Matrice must be squarish for computation of determinant!" );
    
    if( _dimension.columns == 1 )
    {
      return _values[0][0];
    }

    if( _dimension.columns == 2 )
    {
      return (_values[0][0] * _values[1][1]) - (_values[0][1] * _values[1][0]);
    }
    
    double sum = 0.0;
    
    for( int i = 0; i < _dimension.columns; i++ )
    {
      sum += (i % 2 == 0 ? 1 : -1) * _values[0][i] * remove( 0, i).determinant();
    }
    
    return sum;
  }
  
  /**
   * Computes the co-factor of the matrix
   * 
   * @return
   * @throws NoSquareException
   * @throws IndexException
   * @throws MatrixPropertyException 
   */
  public NMatrix cofactor() throws IndexException, MatrixPropertyException
  {
    if( !isSquarish())
      throw new MatrixPropertyException( "Matrice must be squarish for co-factor computation!" );
    
    double d[][] = new double[ _dimension.rows ][ _dimension.columns ];

    for( int i = 0; i < _dimension.rows; i++ )
    {
      for( int j = 0; j < _dimension.columns; j++ )
      {
        d[i][j] = (i % 2 == 0 ? 1 : -1) * (j % 2 == 0 ? 1 : -1) * remove( i, j ).determinant();
      }
    }
    return new NMatrix( d );
  }
  
  /**
   * Computes the inverse of the Matrix
   * 
   * @return
   * @throws NoSquareException
   * @throws IndexException
   * @throws MatrixPropertyException 
   */
  public NMatrix inverse() throws IndexException, MatrixPropertyException
  {
    double factor = 1.0 / determinant();

    return this.cofactor().transposed().multiplyScalar( factor );
  }
  
  /**
   * Returns a boolean value indicating if this matrix is a normal matrix or not
   * 
   * @param m
   * @return
   * @throws IndexException
   * @throws MatrixPropertyException 
   */
  public boolean isNormalMatrix() throws IndexException, MatrixPropertyException
  {
    NMatrix transposed = transposed();
    
    return transposed.multiply( this ).equals( multiply( transposed ));
  }
  
  /**
   * Compare a single row of the matrix with a given row of double values. Both
   * arrays must have the same size, otherwise false is returned.
   * 
   * @param row
   * @param values
   * @return
   */
  public boolean compareRow( int row, double[] values )
  {
    if( values.length != _dimension.rows )
      return false;
    
    for( int i = 0; i < values.length; i++ )
    {
      if( _values[ row ][ i ] != values[ i ] )
        return false;
    }
    
    return true;
  }
  
  /**
   * Provides the values of a particular row
   * 
   * @param row
   * @return
   */
  public double[] getRow( int row ) throws MatrixPropertyException
  {
    if( _values.length != _dimension.rows )
      throw new MatrixPropertyException( "Invalid row index " + row );
    
    double values[] = new double[ _dimension.columns ];
    
    for( int i = 0; i < values.length; i++ )
    {
       values[ i ] = _values[ row ][ i ];
    }
    
    return values;
  }
  
  /**
   * Returns the sum of a row in the matrix
   * 
   * @param row
   * @return
   * @throws MatrixPropertyException
   */
  public double getRowSum( int row ) throws MatrixPropertyException
  {
    if( _values.length != _dimension.rows )
      throw new MatrixPropertyException( "Invalid row index " + row );
    
    double result = 0;
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
       result += _values[ row ][ i ];
    }
    
    return result;
  }
  
  /**
   * Returns a column of the matrix
   * 
   * @param iCol
   * @return
   * @throws IndexException
   */
  public double[] getColumn( int iCol ) throws IndexException
  {
    if( iCol >= _dimension.columns )
      throw new IndexException( "No an valid index, iCol = " + iCol + " matrice does only contain " + _dimension.columns + " columns" );
    
   // double values[] = new double[ _dimension.columns ];
    double values[] = new double[ _dimension.rows ];
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
      values[ i ] = _values[ i ][ iCol ];
    }
    
    return values;
  }
  
  /**
   * Compares an array of double values with a column in the matrix
   * 
   * @param column
   * @param values
   * @return
   */
  public boolean compareColumn( int column, double[] values )
  {
    if( values.length != _dimension.columns )
      return false;
    
    for( int j = 0; j < values.length; j++ )
    {
      if( _values[ j ][ column ] != values[ j ] )
        return false;
    }
    
    return true;
  }
  
  /**
   * Removes a column from the matrix
   * 
   * @param iCol
   * @return
   * @throws IndexException
   * @throws MatrixPropertyException
   */
  public NMatrix removeColumn( int iCol ) throws IndexException, MatrixPropertyException
  {
    if( iCol >= _dimension.columns )
      throw new IndexException( "Invalid column index " + iCol );
    
    double values[][] = new double[ _dimension.rows ][ _dimension.columns - 1 ];
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
      int c = 0;
      
      for( int j = 0; j < _dimension.columns; j++ )
      {  
        if( j == iCol )
          continue;
        
        values[ i ][ c ] = _values[ i ][ j ];
        
        c++;
      }
    }
    return new NMatrix( values );
  }
  
  /**
   * Removes a row from the matrix
   * 
   * @param iRow
   * @return
   * @throws IndexException
   * @throws MatrixPropertyException
   */
  public NMatrix removeRow( int iRow ) throws IndexException, MatrixPropertyException
  {
    if( iRow >= _dimension.rows )
      throw new IndexException( "Invalid row index " + iRow );
    
    double values[][] = new double[ _dimension.rows -1 ][ _dimension.columns ];
    
    int r = 0;
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
      if( i == iRow )
        continue;
      
      for( int j = 0; j < _dimension.columns; j++ )
      {  
        values[ r ][ j ] = _values[ i ][ j ];
      }
      
      r++;
    }
    
    return new NMatrix( values );
  }
  
  /**
   * Provides the upper triangle of the matrix. All elements in
   * the lower triangle are whiped out by setting the value to 0x00
   * 
   * @return Object of type NMatrix 
   * @throws MatrixPropertyException
   */
  public NMatrix getUpperTriangle() throws MatrixPropertyException
  {
    double values[][] = new double[ _dimension.rows ][ _dimension.columns ];
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
      for( int j = 0; j < _dimension.columns; j++ )
      {
        if( j < i )
          values[ i ][ j ] = 0;
        else
          values[ i ][ j ] = _values[ i ][ j ];
      }
    }
    
    return new NMatrix( values );
  }
  
  /**
   * Returns the lower triangle of the matrix
   * 
   * @return
   * @throws MatrixPropertyException
   */
  public NMatrix getLowerTriangle() throws MatrixPropertyException
  {
    double values[][] = new double[ _dimension.rows ][ _dimension.columns ];
    
    for( int i = 0; i < _dimension.rows; i++ )
    {
      for( int j = 0; j < _dimension.columns; j++ )
      {
        if( j > i )
          values[ i ][ j ] = 0;
        else
          values[ i ][ j ] = _values[ i ][ j ];
      }
    }
    
    return new NMatrix( values );
  }
}
